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Abstract. The particle transport of impurities in magnetically confined plasmas 
under some conditions does not find, neither quantitatively nor qualitatively, a 
satisfactory theory-based explanation. This compromise the successful realization of 
thermo-nuclear fusion for energy production since its accumulation is known to be one 
of the causes that leads to the plasma breakdown. 

In standard reactor-relevant conditions this accumulation is in most stellarators 
intrinsic to the lack of toroidal symmetry, that leads to the neoclassical electric field to 
point radially inwards. This statement, that the standard theory allows to formulate, 
has been contradicted by some experiments that showed weaker or no accumulation 
under such conditions BE]. 

The charge state of the impurities makes its transport more sensitive to the electric 
fields. Thus, the short length scale turbulent electrostatic potential or its long wave¬ 
length variation on the flux surface <l>i - that the standard neoclassical approach 
usually neglects - might possibly shed some light on the experimental findings. In the 
present work the focus is put on the second of the two, and investigate its influence of 
the radial transport of We show that in LHD it is strongly modified by 4>i, both 
resulting in mitigated/enhanced accumulation at internal/external radial positions; for 
Wendelstein 7-X, on the contrary, 4>i is expected to be considerably smaller and the 
transport of not affected up to an appreciable extent; and in TJ-II the potential 
shows a moderate impact despite of the large amplitude of 4>i for the parameters 
considered. 


1. Introduction 

From a passive viewpoint the presence of impurities in thermo-nuclear fusion plasmas 
is intrinsic and is bound to the main wall material, whose choice in turn depends on 
the desired wall recycling properties. Furthermore the deliberate introduction of certain 
impurity species is widely recognized as a powerful technique for the reduction of the 
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heat exhaust on the plasma facing components to an tolerable level (see e.g. |3] and ref¬ 
erence therein). The high radiative efficiency of impurities, exploited for that purpose, 
results at the same time in a disadvantage if they accumulate in the plasma core - that 
apart from diluting the fusion fuel - in some cases can lead to the radiative collapse 
of the plasma. The balance between the required concentration in the edge region and 
the lack of its accumulation in the core undoubtedly needs of the identification of the 
physical mechanisms that determine the radial particle flux of impurities, leading them 
to accumulate or to avoid or mitigate that accumulation. 

To that respect, in stellarators and heliotrons the standard physical picture can be sum¬ 
marized as follows. The lack of toroidal symmetry of the magnetic field structure in 
these devices leads to the existence of a population of particles trapped in the helical 
wells that determines the fluxes. Unless that the magnetic field is such that the net 
radial drift of these particles is sufficiently reduced, the total flux surface average radial 
current does not vanish. In standard conditions where the bulk ions and electrons with 
masses and mg have similar temperature Tj ~ Tg their particle diffusion coefficients 
are such that Di ^ ^/mi/rn^Dg jl]. The system needs then of a radial electric field, 
so-called ambipolar, that reduces the ion radial particle flux to the order of that for the 
electrons. This reduction is driven by the confining role that the E x B poloidal preces¬ 
sion related to it has on the bulk ions. Tipically the electron confinement is not affected 
to the same extent due to the more frequent collisions they undergo {ug r\j y^mi/nigi'i), 
which prevent them from completing their poloidal precession orbits. 

Writing the ambipolar electric field as = E^Vr, with Er = —d$o/dr and $o = 'ho(^) 
the electrostatic potential depending by definition only on the flux surface label r, it is 
generally the case that under the previous conditions Er < 0 {ion root regime), or in 
other words that the ambipolar electric field points radially inwards. The importance 
of this for the accumulation of impurities is clear on the context of the standard neo¬ 
classical theory. 

The particle flux density across the flux surfaces for the species a can be expressed 
in terms of the neoclassical transport coefficients Lfj and the thermodynamic forces as 
follows, see e.g. ref. |5], 



( 1 ) 


In eq. (ID. a is the species index, n is the particle density, T the temperature, q and Z 
are respectively the charge and charge state of the species, e is the unit charge modulus 
and the prime represents the derivative respect to r. Therefore it is straightforward to 


conclude that Er < 0 drives particles towards the center of the plasma and out of it 
otherwise. The weight of the charge of the species Z leads inevitably to the dominance 


of the term driving transport through Er as the Z increases, and subsequently supports 
a stronger accumulation of impurities as the ionization state of these is higher. In a 
tokamak the accumulation can be moderated by the so-called temperature screening. 
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This mechanism rests on the fact that the term preceding the temperature gradient in 
eq. can be negative for ions in the long mean free path regime [Imfp). Subsequently 
for peaked temperature prohles the inwards flux of bulk ions driven by T' enforces the 
outwards flux of impurities in order to preserve quasi-neutrality along the radial direc¬ 
tion. Nevertheless this mechanism cannot take place in stellarators since the asymptotic 
scaling of the transport coeflicients with the collisionality, inversely proportional to u 
and ^/l' in the Imfp, is such that Lf 2 lLfi — 3/2 > 0 in both cases. Thus the standard 
prediction for the radial particle transport of impurities in non-axisymmetric devices is 
roughly speaking of accumulation in ion root regime and absence of it in electron root 
plasmas [Hj where Er > 0. 

Nevertheless, the standard neoclassical picture considering only the ambipolar elec¬ 
trostatic potential $0 as constant on each flux surface is particularly questionable in 
non-axisymmetric systems if impurities of moderate to high charge are considered. This 
is due to two reasons: 1) regarding the approximation of the full electrostatic potential 
to only the ambipolar part $0 dismisses a potential variation on the flux surface $1 
that in non-axisymmetric system can be particularly large. This is due to the different 
dependency of fi (first order departure of the distribution function from its equilibrium 
/o) with the collision frequency than that given in axi-symmetric system, and subse¬ 
quently the perturbed linear density ni and related charge density and potential; 2) the 
fact that the weight of the E x B drift and streaming acceleration related to $1 respect 
to the magnetic drift and mirror acceleration scale proportionally to the charge state Z. 

Regarding the first of the aforementioned causes, in an axi-symmetric tokamak the effect 
of the radial electric held is weak and the E^ x B drift is negligible at zeroth order in 
the normalized Larmor radius to the system size drift kinetic expansion parameter: 
p* = pjL. Considering a Fourier-Legendre expansion of fi - Fourier in the poloidal 
coordinate 9 and Legendre in the pitch angle variable ^ = v\\/v - the first order drift 
kinetic equation (DKE) shows that the structure of fi is such that the only non-zero 
Legendre coefficients are the even ones for the cos{m9) terms and the odd for the sin(m6*) 
terms. The latter results in addition proportional to the collision frequency. Thus the 
only terms that can drive a density, and subsequently potential variation on the flux 
surface formally tend to zero in the limit of vanishing collisionality and will have a 
sm{m9) components only. Here m is the Fourier poloidal mode number, ny is the 
parallel velocity and v is the velocity modulus. 

On the other hand the requirement of keeping the Er x B drift at lowest order 
particle trajectory in non-axisymmetric systems breaks that structure in fi. In this 
case the bounce averaged solution of the DKE (/i)^ results in a cos 9 structure with 
amplitude scaling with va/flE- With fhe modulus of the magnetic drift velocity and 
Qe = Er/rBo the poloidal Er x B precession frequency. This picture where a finite 
non-vanishing variation is expected finds a cartoon in the fact that in non-axisymmetric 
system the helically trapped particles shift their poloidal precession orbits from the birth 
flux surface due to the action of both the magnetic drift and the poloidal Er x B 
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Veo 


Vd 


Vd>o X B 


mvj + pB 
~q 


b X V5, 


( 2 ) 

(3) 


with B the magnetic field vector, B its modulus and b = B/5, uy and v±_ the parallel 
and perpendicular components of the velocity and p = v‘j_/2B the magnetic moment. 


The second of the reasons aforementioned about the possible impact of d>i on impurities, 
follows from the comparison between $1 and B as sonrces of transport and trapping. 
Considering the E x B drift related to $1 


vei = 


V$i X B 

52 


it is in comparison to Vd is of the following order, 


(4) 


vei Ze^i R , , 

^ “ T “a’ ^ ^ 

where it has been assnmed that $1 and the toroidal component of B have typical 
variation length-scales VB ~ R~^B and V$i ~ The ratio shown in eq. [5] 

indicates that even for low values of the ratio ed>i/T, $1 can be a source of impurities 
transport of the same magnitnde than the magnetic field gradient and curvatnre. 
Moreover, the streaming acceleration related to $ 1 , Og = —Zeh ■ Vd>i and the magnetic 
mirror am = —ph ■ VB inside the helical wells are of order 


Og Ze^i B 

Um T AB ’ 

where AB is the amplitude of the helical ripple, and it has been assumed that the vari¬ 
ation length-scale of the helical component of B is of order a. Then, the boundaries 
of the trapping regions considering $1 for the case of impurities can substantially be 
modified by those determined by only the magnetic field structure. 


The impact of $1 on the transport has been stndied in the past both analytically and 
numerically mm- The moderate impact found on the bulk species has motivated the 
neglect of $1 in the standard approach for studying their neoclassical transport, but its 
neglect on impurity transport studies has not been fully justified in non-axisymmetric 
system or studies yet apart from recent works [9j. The arguments presented before and 
expressed in the ratios E] and El motivates this work, where the radial particle transport 
of C®"*" impurities is studied including $ 1 . The calculations including $1 requires relax¬ 
ing the nsnal mono-energetic assnmption, and to that end they have been performed 
with the Monte Carlo code EUTERPE for an LHD like equilibrium, one standard con- 
fignration of the Wendelstein 7-X (W7-X) stellarator and a standard confignration of 
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TJ-II. The main features of the code are described in sec. [2] and the results are shown 


in sec. 131 Finally a summary on the amplitude of given in these three devices and 


the conclusions are presented respectively in sections HlEl 

2. The code EUTERPE and the calculation of $1 

The calculations of the radial fluxes for the impurities and $1 has been carried out 
using the Monte Carlo fi particle in cell (PIC) code EUTERPE [TO] . The cnrrent ver¬ 
sion in the gyro-kinetic modality is electro-magnetic, non-linear and considers the full 
surface and radial domain. In the neoclassical version applied to the present problem 
the transport is assumed radially local instead. It performs the collision of the 5f (that 
represents fi here) off the equilibrium distribution function /o, applying at every time 
step a random change in the pitch angle of each marker after the integration of the 
collisionless trajectory The collision frequency for the colliding species a is 

set as the sum of the deflection collision freqnency of a over all the target species b, 
including a itself: Pa = '^ab ^ab- The calculation of under the local ansatz is bound 
to be carried out by iterative adjustment of the valne of Ej. nntil the ambipolarity of 
the fluxes is fulfilled. Since the computational cost is related to the time the fluxes take 
until their convergence and the number of iterations, for the calculations presented in 
sec. 131 Er is a precalcnlated input obtained with a transport code [13] based on DKES 
[I31II3] or with the neoclassical code GSRAKE ng. 

The local ansatz considers that the drifts across the flux surfaces - in this case the 
magnetic drifts and E x B drift vei related to 'hi - are in modnlus mnch smaller - of 
order p* - respect to the parallel velocity, that in tnrn is of order of the thermal velocity 
Uth- The same assumption is considered for the /i part of the distribution function 
respect to the equilibrium /q. This justifies the neglect of the second order terms in 
the left hand side of the kinetic equation • V/i and vei ■ V/i. Since the only drift 
retained at lowest order is the E x B drift related to the ambipolar electric field V£;o, 
and this pushes the markers within the flux surface and not across it, each flux surface 
can be loaded with markers and considered separately of the other integrating in this 
local limit the following characteristic trajectories in the phase space (R, U||,/r): 


m = -b • Vd>i - pb • VR - ^ (b X VR) • V$o, 

m B^ 

/i = 0. 



(7) 


( 8 ) 

(9) 


followed, as aforementioned, by a random kick on the pitch angle of each marker that 
performs the collisional step. Then the Vlasov equation integrated has the following 
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Figure 1: EUTERPE/DKES comparison for the radial particle flux density of C®"'' as 
a function of the effective radius normalized to the minor radius r/a. The plasma 
parameters correspond to the case B.III shown in fig. [3l 


form, 


dt 


+ R - V/i = “/m (vd + V£;i) • Vr 


n 


n 


+ + 


where /m = [no /exp 


— \Vh-\- v\ 


/2u, 


3 Q ^ \ T' 
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2T 2 T J T 

(lb) 

is the standard local Maxwellian. 


To obtain eq. [THl the equilibrium distribution function includes the Boltzmann response 
/o = /Mexp(-Zed>i/r). 

The set of eqs. ITlfTOl in the case $1 were neglected would lead to that of the local and 
mono-energetic limit of the neoclassical theory. In fig. [T] it is shown a radial profile of 
the particle flux density of is such limit, comparing the result obtained using this 
version of EUTERPE with that obtained with a transport code based on the DKES. 
The case corresponds to the LHD-like magnetic configuration considered in the forth¬ 
coming section ITTI for the case labeled as B.III (see figs. EJ, 


Regarding the calculation of $ 1 , instead of the less restrictive ambipolariy constraint 
ZaC (Pq; • Vr) = 0 that determines the ambipolar part of the electrostatic potential 
•ho, the calculation of •hi requires the fulfillment of quasi-neutrality among all the species 
a, 
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= ( 11 ) 

a 

Taking the density for the species a up to first order as the sum of the equilibrium part 
Hao (now varying on the flux surface due to the adiabatic response) and the first order 
departure 


ria = u„o exp(-Z„ed>i/r„) + n„i, (12) 

yields to the following relation in the approximation Ze^i/T <C 1 

T / rp T \ ~^ 

= — ( Uoe + u-oj— + Z'^nQZTpr ) (^n + ^le + Zniz). (13) 

e V Ti Tz J 

For the calculations carried out in this work the electron response has been assumed 
adiabatic. In addition the concentration is set sufiiciently low in all cases (Zeg = 1.1) to 
justify the tracer impurity limit taken into account and neglect the contribution of niz to 
$ 1 . The relaxation of this assumption indeed may conflict with the truncation applied 
to the Taylor series of the adiabatic part so to keep it linear in d>i for the impurity 
distribution function since, as it is shown in next section Z^e^i/Ta can become of order 
unity. Considering more terms in the expansion of the Boltzmann response would lead 
then to a non linear quasi-neutrality equation, whose complexity if left out of the scope 
of this paper. Thust, under the assumptions aforementioned the resulting equation for 
•Fi reduces to. 


$1 


Te f Te 

- I Uoe + Uoi — 

e V Ti 


-1 




(14) 


and reflects that the potential required to balance the lack of charge quasi-neutrality at 
zeroth order is guaranteed by the first order depature from the equilibrium density of 
the bulk ion species, H~^ throughout all the work. 

On the numerical side, the Fourier solver implemented in EUTERPE applies a low 
mode number filter to solve the quasi-neutrality equation that selects the modes in the 
intervals [—4, 4] in the poloidal mode number m, and [0,4] in the toroidal one n. Note 
that $i(6*,0) is obtained at each flux surface separately. As it is different at each of 
them a radial dependency that is not taken into account in the Vlasov equation [10] 
is present. Furthermore, since d>i is obtained at every time step makes that formally 
the dependence of •Fi is more correctly expressed as 'Fi(r, 6*, 0, t). Nevertheless, on the 
one hand the radial dependency can be neglected since the radial electric field related 
to $1 is negligible respect to the ambipolar part: d^i/dt <C d'Fo/dr. And regarding 
the time dependency, this is averaged out in the time interval when the stationary 
conditions have been reached. Then, for the remeinder of this paper we will assume 
the dependence of as it is actually felt by the tracer impurities in our calculations of 
their fluxes, just a stationary potential d>i(d, 0) plugged into eq. dOlto which impurities 
do not contribute given its low concentration. Thus the nonlinear feedback of the 
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Figure 2: (Left) Stationary solution of $1 obtained with EUTERPE as a function of 
the poloidal and toroidal coordinates 9 and <p. Real and imaginary parts (left and 
right respectively) of the Eourier coefficients 'Fpmn on the filtered modes window, and 
normalized to the modulns of the largest component, in this case cos 6. 


Magnetic configurations 


Device 


5o,o(r/a = 0.5) [T] R 


m 


a 


a 


LHD 

1.54 

3.6577 

0.5909 

Wendelstein 7-X 

2.78 

5.5118 

0.5129 

TJ-H 

0.996 

1.5041 

0.1926 


Table 1: Major and minor radii, R and a respectively, and amplitude of the Boozer 
harmonic {m,n) = (0,0) of the magnetic field at the mid plasma radius, Bq q. 


impurity density inhomogeneity on themselves through its impact on $1 is kept out 
of the scope of this work as well. An example of snch map with the corresponding 
spectrnm represented with the absolnte values of the real and imaginary parts of the 
complex Fourier coefficients normalized to the component with maximum amplitude is 
shown in fig. El We advance that these exampie figure correspond to the case B.III 
discussed later in section EH 

3 . Radial flux of C®"*" including $1 

In this section the numerical results are presented and discussed. Three different non- 
axisymmetric devices have been considered: the heliotron type Large Helical Device 
(LHD, Toki, Japan); the helias type Wendelstein 7-X (W7-X, Greifswald, Germany) 
and the heliac TJ-II (Madrid, Spain). For each of them a vacuum configuration has 
been used and the radiai particle flux of G®"*" calculated, comparing the result when 
$1 is neglected with the result when $1 is taken into account. The plasma parameters 
scanned have resulted in 6 different sets of density and temperature profiles for LHD, 
4 for W7-X and 2 for TJ-H. The paramenters related to the magnetic configuration for 
these three devices are listed in table [H 
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Figure 3: Set of profiles considered for the LHD cases. (Left) Electron density profiles. 
(Center) Electron, Hydrogen and C®"*" temperature profiles. (Right) Magnetic field 
modulus in one LHD field period at the radial position r/o = 0.5. 


3.1. LHD results 

Eor the LHD two different density profiles have been considered. A first one corre¬ 
sponding to a high density scenario and a second one corresponding to a low density 
one. These two profiles are represented in fig. [3] (left) and labeled respectively as A (high 
density) and B (low density). The equilibrium density Uoj of the bulk ions (Hydrogen 
nuclei) and impurities hqz - in all cases fully ionized Carbon C®+ - is determined so to 
fulfill quasi-neutrality among these and the electrons: Z^enoa = 0. The reference 

value considered for the effective charge is for all cases =1.1, that allows to assume 
negligible the impact of the radial flux of Carbon on the ambipolar electrostatic po¬ 
tential $ 0 , and of the influence of the impurity perturbed density niz on the potential 
$ 1 . Three temperature profiles for C®+, set equal to the bulk ion temperature T*, have 
been used, letting the electron temperature Tg profiles fixed. These profiles, sorted by 
increasing ion temperature are labeled as /, II and III, and are shown in fig. |3] (cen¬ 
ter). Finally the magnetic field strength of the magnetic configuration is at the position 
r/a = 0.6 illustrated in fig. |31 (right) as a function of the poloidal and toroidal angular 
coordinates 6 and (f). 

In fig. IHfor the low density set of profiles A the electrostatic potential depending on the 
poloidal and toroidal coordinates d>i(0,(/)) is represented on the top row for the three 
mentioned ion and C®"*" temperature profiles I, II and III. Equivalently the bottom row 
contains for the low density set B the three maps of d>i for the same three temperature 
profiles. All the six maps are considered at the radial position rfa = 0.6. 

It can be observed the qualitative trend mentioned in the introduction. As the maps 
are viewed from left to right and top to bottom, this is with monotonically decreasing 
ion collisionality, the variation of the potential from peak to peak is increased in one 
order of magnitude, accordingly with the decrease of the collisionality (of a factor 40 
approximately). This justifies to some extent the neglect of the electron contribution to 
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Figure 4: At the position r/a = 0.6 for LHD: on the top row $i(0, (p) for the high density 
case (A) and on the bottom row for the low density case (B). The left, center and right 
column corresponds to the temperature profiles considered I, II and III respectively. 


•Fi we are considering, since in the cases considered the electron to ion collisionality is 
also of the same order larger than that of the bulk ions. 

Regarding the central qnestion of how strongly is the impnrity transport affected by •Fi, 
the potential was calculated at the radial positions rfa = {0.2, 0.4,0.6, 0.8} for each set 
of profiles, and the resulting maps of $1 considered as an input for the calculation of the 
C®’*' radial particle finx (F • Vr). The radial particle finx density of as a function 
of the flux surface label is shown in fig. El for each of the 6 profiles sets considered in 
this section, and normalized to the equilibrinm impurity density Uoz- The dotted line 
represents the case where only the ambipolar electric field is taken into account, and 
the solid onee corresponds to the case where both the ambipolar part of the potential 
$0 and $1 are both taken into account. Each of the 6 plots is accordingly displayed at 
the same position than the maps in fig. ID do refered to the pair of profiles considered. 
This is, on the top row the high density scan A and on the bottom one the low density 
one B are shown; and from left to right the figures corresponds to the Ti = Tc(^+ profile 
/ (left), II (center) and III (right). 

At a first glance, reading the plots of fig. |5] from left to right and top to bottom and 
looking at the maps in 01 as orientative for the order of magnitude of $1 for each case 
(note that the maps shown are for the position r/a = 0.6 only), it is found that the 
impact that $1 prodnces on the ratial particle finx of Carbon is continnously increasing 
as d>i does. And interestingly, the modification that this represents respect to the 
standard neoclassical pictnre does not point always to the same direction. <Fi both 
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Figure 5: Radial flux density of C®"*" as a function of r/o inclnding d*! (solid line) and 
neglecting it (dotted line). From left to right the temperature profiles I, II and III are 
considered and from top to bottom the high A and low B density cases. 


mitigates and amplifies the trend of C®’*' to accumulate. On the one hand, in the high 
density scan (three figures on top) the impact of results in a weak mitigation of the 
inwards flux at all radial positions considered, almost negligible in the lowest Tj case 
(top, left) but appreciable in the highest Tj one (top, right). On the other hand, the 
three figures at the bottom corresponding to the low density temperatnre scan coincide 
to the fact that •Fi reduces the inwards Carbon flnx at the internal radial positions 
below r/a = 0.4 approximately and makes it stronger from there ontwards, resnlting in 
the less collisional of all the cases considered B.III (bottom, right) in the reversal of the 
behavior predicted by the standard neoclassical approach and a snbstantial redaction 
of the inwards flux. 

In order to sketch qualitatively what underlies the fact that the radial flux predicted 
by the standard neoclassical prediction becomes stronger or weaker when d*! is taken 
into account, in fig. El the real and imaginary part of the complex Fourier coefficients of 
•Fi, or equivalently the amplitnde of cosine and sine components, is presented at each 
position simulated in the set of figs. El Each of the plots in figs. El refers to that at the 
same position in figs. El All the cases point ont to the fact that the main components 
are the cos 9 and sin 9, and that a correlation between the sign of the latter, plns/minns, 
corresponds respectively to the situation where less/ more inwards flux is driven. This 
is remarkably clear in the low density temperatnre scan. Looking at the three plots 
at the bottom of the figures El and the corresponding Fourier coefficient at the bot¬ 
tom row of figs. El the transition from a mitigated to an enhanced inwards flux accurs 
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approximately at the radial position where the sign of the sin 0 component flips the sign. 


It is important, in oder to provide an idea of the more complex coupling of V$i (in 
comparison with V5) to /i and produce transport, the presence of sine components 
in the spectrum of d*!. This point out to the lack of the stellarator symmetry that for 
instance the magnetic field modulus B has. It is indeed the stellarator symmetry of B 
what makes a necessary condition its absence in $i. Considering the radial flux density 
of particles expressed as 


{r-Vr) = Qd\^ry (15) 

where (...) is the flux surface average. For the bulk ion neoclassical flux it can be written 
as dr/dt ~ Vd • Vr. This term introduces derivatives respect to 9 and (p in the integrand 
of [IS] that makes that a stellarator symmetric magnetic field does not drive any transport 
unless the distribution function indeed develops sine components. This sine components 
are then inherited by the moments of the perturbed part of the distribution function 
like the perturbed density and subsequently the electrostatic potential $i. Regarding 
the impurity radial particle transport one needs to recover the initial hypothesis that 
retaining vei is necessary and express dr/dt ~ (vd + vei) • Vr. Then, regarding how 
$1 participates on the transport of impurities becomes less trivial than for the case of 
B, given that Fourier spectrum has both sine and cosine parts that drives particle 
flux coupled to the perturbed distribution function through both its cosine and sine 
components respectively, and not only through the latter ones as does. To this 
respect the question of how efficient can a specific component of be in counteracting 
the trend of the impurities to accumulate depends on the spatial dependence of their 
distribution function, and brings to the same importance the the knowledge of and 
the impurity density variation on the flux surface. An indeed, both magnitudes have 
recently met an increasing interest in their measurement nziiiHiiin]. is kept out of the 
scope of this paper will be addressed in future works. 

3.2. Wendelstein 7-X results 

The main features of the dependence of with the plasma parameters, as well as the 
role of its Fourier spectrum on supporting or reducing the tendency of to accumulate 

have been discussed in sec. EH In this one we show similar numerical simulations for 
the Wendelstein 7-X stellarator, and a set of profiles expected under operation. These 
are shown in fig. [71 The electron density profile is represented in fig. |7| (left), and as 
before Z^s = 1.1 is considered to determine noj and Hqz. The temperature, the same 
for all species in each of the two cases, is the parameter scanned. Four different profiles 
have been considered. The temperature profiles in fig. Ej (center) are labeled as I, II, III 
and IV increasingly with T. 
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Figure 6: Amplitude of the main Fourier components as a function of the radial position 
for the six cases considered, displayed as in figs. El 
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Figure 7: Set of profiles considered for W7-X. (Left) Electron density profile. (Center) 
Electron, bulk ion and temperature profiles, assumed the same for all the three 
species. (Right) Magnetic field modulus in one W7-X period at the position rja — 0.6. 


The electrostatic potential is in this case appreciably much lower than at LHD. Although 
the collisionality does not reach in any of the the four cases considered in this section 
the value of the lowest collisional low density LHD set R, the variation of the potential 
is noticeably much lower in W7-X for a similar collisionality. 

A possible reason underlying the low electrostatic potential variation in W7-X is the 
physical target for the design of this device, namely, reducing the neoclassical losses and 
bootstrap current by approaching an omnigeneous magnetic field structure. In a perfect 
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Figure 8: At the position r/a = 0.6 for W7-X, $i(6*,0) for the temperature profiles /, 
II,III and IV (from left to right). 


omnigeneous magnetic field the bounced-averaged drift orbits of the localized particles 
align to the flux surfaces. Although this ideal situation is not, and cannot be the 
case, the closer the W7-X magnetic field structure is from it - compared to the other 
devices considered in this work - certainly leads to a smaller departure of the localized 
particles from the flux surfaces and results in a weaker electrostatic potential variation 
on them. Regarding the Fourier spectrum of $i, similar features to LHD’s at similar 
collisionality (cases A.I-III) are shown. A dominant sind component, and a weaker cos 9 
one that increases as collisionalily decreases is the footprint of the electrostatic potential 
in both W7 — X and LHD, although the rest of the modes represent a more noticeable 
contribution to the spectrum of $1 than in LHD, where they are all marginally zero 
with a some exception (e.g. the cos 26*). 

The radial profiles for the particle flux density of C®'*' are represented in fig. [HI in the 
three W7X cases, increasing the temperature from the left to the right. Not surprisingly 
the difference that $1 makes respect to the standard neoclassical prediction can be well 
considered for the present set of parameters and impurity species negligible and aligned 
to the neoclassical optimization aforementioned. 

3.3. TJ-II results 

Finally in the present section two cases are considered for the TJ-II stellarator. Two 
density profiles, shown in fig. fTni lefth at a fixed temperature for the electrons and bulk 
and impurity ions {H~^ and as in the previous cases) represented in fig. flOlf center) . 
These profiles correspond to typical NBI-heated TJ-II plasmas. The contour plot of the 
magnetic field at the flux surface r/a = 0.6 is represented in fig. mUf rightb 
In fig. [H] (top) the electrostatic potential 'Ll is represented for the two cases consid¬ 
ered together with some Fourier harmonics (bottom). The calculation of $1 performed 
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Figure 9: Radial flux density of C®’*' as a function of r/a including (solid line) and 
neglecting it (dotted line). The cases / to IV are represented from left to right. 



Figure 10: Set of profiles considered for TJ-II. (Left) Electron density radial profile. 
(Center) Electron, bulk ion and temperature profiles, assumed all equal. (Right) 
Magnetic field modulus in one TJ-II period at the position rfa = 0.6. 


shows a low $1 amplitude compared to the LHD cases, and of the same order of magni¬ 
tude of a few volts as in W7-X. Nevertheless, the much lower temperature (and higher 
collisionality) in this TJ-II sample compared to those considered for the previous two 
devices, leads indeed to conclude that TJ-II exhibits the strongest potential variation 
of the three, specially if the electrostatic energy variation related to is compared to 
the thermal one. A discussion on this comparison for all the cases considered in the 
present work is given in section HI We advance that the electrostatic potential variation 
for the bulk ions reaches up to a 5-8 % of the thermal energy for the low density TJ-II 
case and up to 10-15 % for the low density one, but despite of this the radial transport 
when $1 is included is not as different from the prediction without it could be expected. 
This can be observed in figs fTW left) and fTW right) where the comparison of radial par¬ 
ticle flux density with and without is represented for the high and low density cases 
respectively. 

A possible explanation could rest on the broad spectrum that 'Ll exhibits in TJ-II, and 
that could lead to a partial cancellation of the effect of Eourier components contributing 
both to enlarge and diminish the flux level. In fig. HD (bottom) the main Eourier 
components are respectively represented as a function of r/a for the high density case 
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Figure 11: At the position r/a = 0.6 for TJ-II standard configuration: d>i(d,(/)) for the 
temperature profiles I (left) and II (right). 



r/a 



r/a 


Figure 12: Radial particle fiux density of C®’*' as a function of r/a including d>i (solid 
line) and neglecting it (dotted line), for the high and low density cases (left and right 
respectively). 


(left) and low density one (right), showing the absence of a clear dominant component. 
Although the cos 9 term is large along the radial coordinate in both cases, its amplitude is 
not much larger than other terms similarly large, as e.g the sind, the helical components 
cos(d + 0), sin(0 + (j)) or sin0. 

4. Final remark on the amplitude of 

In fig. [in] the maximum variation of electrostatic energy for H"*" related to the difference 
between the maximum and minimum value of the potential on the surface APhii is 
represented, as a function of the normalized collision frequency normalized to the bounce 
frequency, u* = r'i,th/<^6, with Ub = Vi\^i/R. This ratio eA^i/T is represented for all 
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Figure 13: Normalized peak-to-peak variation of the potential A^i for all cases 
considered for LHD, W7-X and TJ-II as a function of the normalized thermal collision 
freqnency of bulk the ions normalized to the bounce freqnency coi,- 


the cases considered in the work at the four radial positions considered in each case. 
The value eA'Fi/T provides an estimate of how reliable is the standard neoclassical 
assumption of mono-energetic trajectories for the bulk ions, and if multiplied by the 
corresponding charge state, for the impnrity species of interest. 

To this respect, having taken C®'*' as the impurity for the present study, it is clear that 
the mono-energetic assumption breaks down in the LHD and TJ-II cases at the most 
external position r/a = 0.6 and 0.8, that correspond to the nppermost points of each 
curve. This is remarkably clear at the low collisional LHD set of profiles and the two 
cases in TJ-II. In TJ-II the potential variation is particnlarly large considering the much 
larger collisionality than in the LHD set B, as it was noted at the end of the previous 
section. At similar collisionality, comparing the LHD set case A.I and TJ-II low density 
one, LHD shows roughly an order of magnitude lower 'Ll than TJ-II. Note that in this 
normalized picture the ratios comparing the strength of Vd>i and as transport and 
trapping sources in eqs. [5] and |6] become of order nnit as well. 

In the figure it can be appreciated how at similar collisional regime W7-X reaches values 
of eA$i/r one order of magnitnde lower that in the LHD high density set, and finds 
the usual neoclassical neglect of d>i well valid for the parameters considered. 

5. Summary and conclusions 

The impurity transport is known to be strongly influenced by the ambipolar part 
of the electrostatic potential, nsnally considered as an approximation to the fnll 
neoclassical electrostatic potential. In the present work, we have also considered the 
potential variation on the flnx surface after solving the quasi-neutrality equation, 
and calculated the impact on the radial particle flux of for the devices LHD, W7-X 











$1 and its impact on impurity transport 


18 


and TJ-II, using the code EUTERPE. 

The calculation have confirmed the initial assumptions on the importance of taking 
$1 into account for the transport of impurities, given that for these $1 can both 
modifying substantially the topology of the trapped/passing boundaries or driving 
transport through the related E x B drift in the same degree than VB. 

In most cases $1 has led to appreciable deviations from the results framed in the 
standard neoclassical approach, both mitigating or enhancing the trend to accumulate. 
This was particularly clear in the LHD low collisional cases, and in less extent in the 
higher collisional TJ-II ones, despite the large variation of $1 the latter device exhibit. 
Finally in W7-X $1 has resulted to be weak enough to be negligible for the equilibrium 
and parameters considered. 

Although an extension of previous studies, and apart from a more exhaustive parameter 
and magnetic configuration scan, there are still major points to address and assumptions 
to relax in next works: first of all the relaxation of the tracer impurity limit, since 
scenarios with accumulation unavoidably can lead in long discharges to sufficiently high 
impurity concentration to make this approximation break down; once impurities are no 
longer tracers and must be included in the quasi-neutrality equation the linear adiabatic 
response considered, specially if they have large charge states, is no longer valid either; a 
detailed study of how $1 and the spacial impurity distribution function couples to each 
other could also help to identify what phase shift between both can bring mitigated 
accumulation scenarios. 
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